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A system of N classical particles in a 2D periodic cell interacting via a long-range attractive potential is 
studied numerically and theoretically. For low energy density U a collapsed phase is identified, while in the 
high energy limit the particles are homogeneously distributed. A phase transition from the collapsed to 
the homogeneous state occurs at critical energy Uc- A theoretical analysis within the canonical ensemble 
identifies such a transition as first order. But microcanonical simulations reveal a negative specific heat 
regime near Uc. This suggests that the transition belongs to the universality class previously identified by 
Hertel and Thirring (Ann. of Physics, 63, 520 (1970)) for gravitational lattice gas models. The dynamical 
behaviour of the system is strongly affected by this transition : below Uc anomalous diffusion is observed, 
while for U > Uc the motion of the particles is almost ballistic. In the collapsed phase, finite A''-effects act 
like a "deterministic" noise source of variance 0{1/N), that restores normal diffusion on a time scale that 
diverges with A''. As a consequence, the asymptotic diffusion coefficient will also diverge algebraically with 
N and superdiffusion will be observable at any time in the limit N — > oo. A Lyapunov analysis reveals 
that for U > Uc the maximal exponent A decreases proportionally to N"^^"^ and vanishes in the mean-field 
limit. For sufficiently small energy, in spite of a clear non ergodicity of the system, a common scaling 
law A oc U^^^ is observed for any initial conditions. In the intermediate energy range, where anomalous 
diffusion is observed, a strong intermittency is found. This intermittent behaviour is related to two different 
dynamical mechanisms of chaotization. 
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I. INTRODUCTION 

Theoretical and computational studies ||^,^ have 
been devoted in the last decades to the study of thermo- 
dynamical properties of iV-body gravitational Hamilto- 
nian models. One of the main results is the identification 
of a phase transition from a collapsing phase (CP) ob- 
served at low energy to a homogeneous phase (HP) at 
high energy. In the CP a fraction of the particles form 
a single cluster floating in a homogeneous background of 
almost free particles. Above a critical energy Uc the clus- 
ter disappears and all the particles move almost freely. In 
the transition region, the specific heat becomes negative 
in the microcanonical ensemble. This instability has been 
widely studied in astrophysics, where it is known as the 
gravothermal catastrophe [^-^. A negative specific heat 
regime seems to be thermodynamically inconsistent. But 
this paradox is solved once the non-equivalence of canon- 
ical and microcanonical ensembles in the neighbourhood 
of such phase transition is demonstrated ||l|,0|. These re- 
sults have been sucessfully confirmed by numerical inves- 
tigations of self-gravitating non-singular 2D systems with 
short range interaction [pP] and by the study of a sys- 
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tern of particles confined on the surface of a sphere with 
variable radius The reported analysis were usually 
focused on 2D systems, but the features of this thermo- 
dynamical transition are expected to be independent of 
the space dimensionality |jl[|. 

More recently it has been shown that cluster formation 
in simple ID iV-body models with long-range attractive 
interaction |8|J^ is similar to the Jeans instability oc- 
curring in gravitational systems pd| . In particular, fol- 
lowing this analogy Inagaki and Konishi [ pdl were able 
to derive the correct critical temperature for the declus- 
tering transition ]l^ . 

Preliminary indications suggesting that clustering phe- 
nomena can have effects on single particle dynamics have 
been reported for 

• (i) a system of fully and nearest-neighbour coupled 
symplectic maps with attractive interaction |13| ; 

• (ii) atomic clusters Q; 

• (iii) turbulent vortices [p^ , p^ . 

All these systems exhibit a clustered phase associated to 
an anomalous diffusion law. The anomalous transport 
is related to the single particle dynamics that intermit- 
tently shows a sequence of localized and almost ballistic 
behaviours, corresponding to the trapping and untrap- 
ping in the cluster, respectively. 

Anomalous diffusion has been revealed in dissipative 
and Hamiltonian models |17 1^ as well as in experi- 
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merits The first investigated example of anoma- 

lous behaving systems was a one dimensional chaotic 
map proposed by T. Geisel and co-workers JTof that was 
introduced to describe enhanced diffusion in Josephson 
junctions poj . More recently, anomalous diffusion was 
also studied in the context of fluid dynamics and solid 
state physics both experimentally |^,^,^ and theo- 
retically ||l^,|l8|,0. However, the conservative models 
usually studied concern the motion of a single particle 
that evolves in a fixed two-dimensional periodic potential 
p2| , p7 23 2^ or the dynamics of low dimensional sym- 
plectic maps. Only few studies have been devoted to 
extended models with iV >> 1 p|,p^j28|. 



The main purpose addressed in this paper is to estab- 
lish a link between equilibrium results, transport proper- 
ties and finite A^-effects for A^-body Hamiltonian models 
with long-range attractive interaction. In order to ex- 
amine in detail these points we consider a simplified 2D 
self-gravitating model, where the dynamics of each parti- 
cle can be written in terms of mean-field quantities that 
are self-consistently defined through the coordinates of all 
the particles. Therefore each particle is moving in a non 
autonomous self-consistent potential depending on mean- 
field quantities. This peculiarity allows to implement a 
very efficient code and therefore to examine systems with 
a high number of particles {N ~ 4, 000 — 10, 000) a for 
long time periods {t ^ Ix 10^ - 1 x 10^, where the typical 
time of evolution is 0(1)). 

We observe a declustering transition, analogous to 
the one studied by Hertel and Thirring, and, in corre- 
spondence of the CP, anomalous diffusion is detected. 
Above the transition point also the dynamical behaviour 
changes and the particles move almost freely. 

In the CP, the shape of our self-consistent single par- 
ticle potential resembles that of the so called egg-crate 
potential. Several studies were devoted to the motion of 
a single particle in an autonomous egg-crate potential, 
i.e. when the potential represents a fixed landscape for 
the particle 0,|8 22-2^. In particular, anomalous dif- 
fusion was observed and such behaviour was explained 
as due to trapping and untrapping of the particle orbit 
in a self-similar hierarchy of cylindrical cantori As 
a result the trajectory of the particle on the 2D surface 
is similar to a Levy walk. This suggests that anomalous 
transport, present in our model, should be explained by 
similar mechanisms. 

Due to the self-consistent nature of the potential in the 
model we consider, finite A^-effects induce a time fluctu- 
ation in the potential seen by a single particle and essen- 
tially act as a noise source of typical intensity 0(1 /^/N). 
Environmental fluctuations are known to restore normal 
diffusion on long time scales when added to a dynami- 
cal systems exhibiting anomalous diffusion |^,^. In- 
deed, this is what happens also in our model, apart that 
now the fluctuations are intrinsically related to its deter- 
ministic dynamics and not to an external bath. In the 
mean-field limit the fiuctuations disappear and anoma- 
lous diffusion is present at any time. 



The study of the degree of chaoticity in the present 
model reveals that the system (for low energy) is highly 
non-ergodic and, despite the fact that thermodynamical 
properties seem not to depend on initial conditions, the 
dynamical indicators (e.g. the maximal Lyapunov expo- 
nent) are heavily affected from the initial state of the sys- 
tem. Despite this non-ergodicity, the low energy depen- 
dence of the maximal Lyapunov shows a common scaling 
law for all the considered initial conditions. In the inter- 
mediate energy range, two distinct chaotic mechanisms 
are identified in our model and give rise to intermittency 
in the dynamical evolution. For high energy the system 
becomes integrable in the mean-field limit and the Lya- 
punov vanishes as an inverse power law of A^. 

This article is organized in the following way. In sec- 
tion II, we present the Hamiltonian model we focus on. 
In section HI, we report a thermodynamical description 
of the model within the canonical ensemble. Equilibrium 
values of the mean-field are derived and the existence of a 
phase transition from a homogeneous to a collapsed phase 
is shown. Section IV is devoted to numerical description 
of dynamical properties. In particular, we'll concentrate 
on the diffusion of single particles and on the inffuence 
of finite A^-effects. The chaotic behaviour of the model is 
discussed in section V in a broad energy range. Finally, 
in section VI the reported results are briefiy discussed. 



II. THE MODEL 

In the present paper we study the static and dynamical 
properties of a A^-body system enclosed in a 2-D periodic 
cell, with a classical long-range inter-particle potential. 
The dynamics of each particle is ruled by the Hamiltonian 
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3 — cos(a;i — Xj) 



' cos(2/i - %) - cos(a;i - x^) cos(?/i - yj) 



= K + V 



where (xi^yi) G [-tt, tt] x [-7r,7r], {xi,Pi^x) and {yt,Pi,y) 
are the two couples of conjugate variables. The particles 
are assumed to be identical and to have unitary mass. 
K (resp. V) is the kinetic (resp. potential) energy. The 
reference energy is chosen in such a way that the energy of 
the system vanishes when all the particles have the same 
position {V = 0) and zero velocity {K = 0). The presence 
of the third term in the potential energy is essential in 
order to ensure an interaction between the two couples 
of conjugate variables. 

The considered interparticle potential belongs to the 
following class of 2— D periodic potentials : 
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3 = 1 0<k?<s 



^ c(fc)[l -cos[k« (r-rj)]j (2) 



where s is a parameter that determines the number 
of harmonics included in the Fourier expansion of V , 
r = {x,y) and k — {n,j;,ny) is the wave-vector, with 
rix and ny two integers. c(fc) = c{~k) is a real valued 
function of the modulus of k, that fixes the coupling 
strength of the harmonic k. For model (]l|), c{k) = 
and s = 2. For the same choice of c(fc), but in the limit 
s oo, a 2 — D self- gravitating Newtonian potential 
Voo{r) cx X]t=i l°Sk ~ is recovered, once a rescaled 



time t/y/N is considered. It should also be noticed that 
interactions of the type log \ r\ among particles arises also 
in point vortices model for 2D turbulence There- 
fore the simplified model ([l]) should share some common 
behaviour with point vortices systems. 

The interparticle potential appearing in (||) is essen- 
tially the Fourier expansion of Voo (r) limited to its first 
three terms and the corresponding equations of motion 
for the coordinates 



1 ^ 

N ^ 

j=i 



iin{x, - Xj) + sin(xi - Xj) cos{y, - yj) 



(3) 

due to the symmetry of (0), the equations of motion for 
y are obtained exchanging x <-> y. 

Previous investigation of model (|l|) |^ has revealed 
that at low energy U = H/N the particles, due to at- 
tractive coupling, are organized in a unique collective 
structure {cluster) with {/-dependent spatial extension. 
This clustered phase (CP) survives up to a critical en- 
ergy Uc ~ 2. Above Uc a phase transition occurs to 
a homogeneous phase (HP). This transition will be dis- 
cussed in more detail in the next section. Here, we will 
limit to evidence collective behaviours in the system by 
rewriting model (^ in terms of the following mean-field 
variables : 



Mi= (< cos(a;) >Ar, < sin(a;) >m) 
= Mi (cos($i),sin($i)) 



(4) 



M2= (< cos{y) >N, < sm{y) >n) 

= M2 (cos($2),sin($2)) (5) 

Pi= (< cos(a; + y) >n,< sm{x + y) >jv) 

= Pi (cos(*i),sin(*i)) (6) 

P2= (< cos(a; - y) >n, < sin{x - y) >n) 

= P2 (cos(*2),sin(*2)) (7) 

where < .. >7v denotes the average over all the particles 
in the system; (j)z (resp. Vz) and (resp. Pz) are the 



argument and the modulus of the mean-field vectors 
(resp. P2) with z = 1,2. The moduli AI^ and Pz are 
maximal and equal to 1 when all the particles have the 
same position and their value decreases when the spatial 
distribution of the particles extends. For an homoge- 
neous distribution, due to finite N effects and according 
to the central limit theorem, Mz Pz k. 0(l/y/N) for 
z = 1, 2. Therefore, the quantities Mz and Pz can be 
thought as order parameters characterizing the degree of 
clustering of the system ]l3|,pt. 

By reexpressing the equation of motion (|^) for both 
coordinates x and y in terms of the mean-field variables 
Mz and Pz, one straighforwardly shows that the time 
evolution of each particle i is ruled by the following single 
particle non-autonomous Hamiltonian : 



'^Pi cos{xi +yi~ iIji) - ip2 cos(2:i -yi- 1P2) 



3 — Ml cos{xi — — M2 cos(?/i — ^2) 

(8) 



where Ki and Vi are the 1-particle kinetic and poten- 
tial energy, respectively. Hence, each particle moves in a 
mean field potential Vi determined by the instantaneous 
positions of all the other particles of the system. 

Since V is invariant under the transformations x <-> 
—x, y <-> —y and x <->?/, it turns out that in the 
mean- field limit {N 00 with constant U = H/N), 
Ml = M2 = M and Pi P2 = P- Moreover, in this 
limit and assuming 4>z = ipz = 0, the 1— particle potential 
is similar to the egg-crate potential studied in references 
2^ . In Fig. (|^) the shape of Vi at a given time 
1.00 and iV = 4, 000 is shown. The potential is 
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for U 



periodic along the two spatial directions and in each ele- 
mentary cell it has 4 maxima (Vm = 3-I-2M— P), 4 saddle 
points {Vs = 3 + P) and a minimum {Vm = 3 — 2M — P). 
The depth of the potential well is ~ {Vs — Vm) and the 
center of the cluster will coincide with the position of the 
potential's minimum. In the limit t/ — > 0, M and P — > 1 
and therefore Vm ^ 4, 14 ^ 3 and Vm — > 0. In this limit- 
ing situation all the particles are trapped in the potential 
well they create. For increasing U, the kinetic contribu- 
tion becomes more relevant and the number of clustered 
particles reduces. This implies that also the values of M, 
P and of the well depth decrease. For U > Uc, the par- 
ticles are no more clustered and in the limit N ^ 00 the 
single particle potential becomes flat Vm — Vs ^ Vn ~ 3 
and time independent. 

However, due to finite N effects, the instantaneous 
mean-field variables fluctuate with typical time 1 within 
the statistical band range ~ 1/\/N. This implies that in 
the CP the single particle will move in a fluctuating po- 
tential landscape. As a consequence, the system admits 
a time pulsating separatrix that sweeps a phase space do- 
main proportional to 1/ in the neighborhood of the 
average position of the separatrix. 

Due to the shape of V^, we can qualitatively distinguish 
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at each time three sub-groups of particles: the low energy 
trapped particles (LEP) oscillating in the self-consistent 
well of the potential (with energy hi Vm), the interme- 
diate energy particles (lEP) evolving in the stochastic sea 
neighbouring the separatrix (hi ~ Vs) and the high en- 
ergy particles (HEP) moving almost freely (hi >> Vm)- 
Obviously, during the time evolution each particle can 
pass from one sub-group to another. 



III. THERMODYNAMICAL ASPECTS 

A. Canonical Prediction for the Equilibrium 
Properties 

The problem of finding the thermodynamical poten- 
tials of model (|]) in the canonical ensemble is resolved 
starting as usual from the partition function Z{N,f3), 
where /3 is the inverse temperature. For Hamiltonian (|^), 
the partition function factorizes in Z{N,P) x Zy, 

where : 



where G = {u 



:,)/ f3 — logi? -I- log/3 is a strictly pos- 
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Zy OC 



exp 



f3— (2M2 + 



exp(cyl^) cx / exp(— -I- 2y/cuA ■ A)(i^UA 



(iN/2) and perform the rescaling ^/2P/NUl„ — > u^, 
w = M, P. Then expressing u^, — [uw cos 9^ , Uw sin ( 
in polar coordinates, we obtain: 
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is the kinetic part obtained assuming a Maxwell- 
Boltzmann distribution for the momenta. The potential 
contribution Zy is evaluated by considering the expres- 
sion of the potential V in terms of mean-field variables 
and P^. As already mentioned, in the mean-field 
limit due to the symmetries of V we have Mi = M2 = M 
and Pi = P2 = P- Consequently, V can be reexpressed 
as a quadratic form in M and P and the resulting ex- 
pression for Zy reads jo) : 



d^x d^y (10) 



where S stands for the unit circle and where the constant 
part of V is omitted. We can evaluate the integral in Eq. 
( p^ ) by using the Hubbard-Stratonovich trick: 



(11) 



where A is a two component real vector and c a positive 
constant. Following PJs^], we substitute expression ( [ll| ) 
into Eq. (|l^ with {A = M,c = l3N) and {A = P,c ^ 

for 



itive real valued function independent of N and R. The 
main contribution to Z in the limit N ^ 00 can be 
evaluated with the saddle point method, and reduces to 
the estimation of the minimum of the function G. The 
minimum value will not depend on 9m and 9p, because 
R is maximal with respect to these two variables for 
{9t,[,9p) = (0,0) for any um and up. Its expression 
thus reduces to: 

Rm{uM,up) = max R{um,up,9m,9p) = 

Om ,9p 

= 27r y Io(^um + V2upcos{s)'^ cxp^UAi cos{s)^ds (13) 

where Iq is the zero order modified Bessel function. Fi- 
nally, the minimum of G = G(/3, um, up) is obtained nu- 
merically using expression (|l^) and a Raphson-Newton 
scheme. The (3 dependence of the solutions {uM,up) for 
which G is minimal is reported in Fig. |[ At large (3 
(low temperature) the minimum is located far away from 
the origin and corresponds to the CP i.e. to the broken 
symmetry phase. When f3 — 2, a saddle-node bifurcation 
occurs at the origin {um,up) = (0,0) (corresponding to 
the HP), where a second minimum appears. This can be 
shown with a Taylor expansion of G in the neighborhood 
of the origin : 



(14) 



From the above equation, one can observe that for /? > 2 
the origin is a saddle point (it is a local maximum if 
/3 > 4) while for /3 < 2 the origin is a local minimum. The 
upper inset of Fig. |2| shows the details of the temperature 
domain where both minima coexist. When decreasing /3, 
a second bifurcation occurs for (3 = 1.81 and the mini- 
mum corresponding to the broken symmetry disappears 
while the one at the origin persists. 

The mean-field variables (^,(^),(|6|) and (^ can be es- 
timated from their joint characteristic function: 



G{P,UM,Up)Ki Um 75-0 



-log Uvr^ +ln/3 
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*(fTM,Tp) =< e('^M'M+ap.P) 



> 



(15) 



where the average < • > runs over the corresponding 
Gibbs ensemble. The quantities ctm and crp can be con- 
sidered as external fields while M and P are the equiva- 
lent of magnetizations or polarization vectors. The func- 
tion ( p!5| ) can be explicitly computed and its first mo- 
menta are the moduli of the mean-field variables : 



M= 



2n 



d\aM\ Rm{uMiUp) 

/ UMUpexp{-NG) duMd9Mdupd9p (12) / cos(s)/o (um + \/2up cos(s) j exp ma/ cos(s) (16) 



Z = Zk X Zy CX 
'(KxS) 
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p= 



27r 



d\(Tp\ Rra{uM,Up) 

ds cos(s)/i {um + "\/2mp cos(s)^ exp um cos(s) (17) 

where function Rm was introduced in Eq. (^), Ii is the 
first order modified Bessel function. From these expres- 
sions it is easy to check that when {um, up) coincide with 
the origin, then M = P ~ 0. This aUows us to identify 
such solution with the HP. While if {um,up) is located 
away from (0,0) the solution corresponds to the CP. 

More specifically, M and P can be considered as two 
order parameters that describe the spatial extention of 
the cluster in the {x, y) plane (i.e. the degree of clustering 
of the system). The analytical estimation for M and P 
are plotted in Fig. |3|(a) (solid line) as functions of U . 
For C/ ^ 0, M and P tend to 1 as expected for a fully 
ordered phase. For increasing C/, kinetic effects increase, 
the size of the cluster thus grows and consequently M and 
P decrease. For l.Q <U < 2.0 the analytical values are 
not reported, because in this energy range it is impossible 
to estimate them with a canonical approach. This simply 
reflects the first order nature of the transition. For U > 
2.0, M = P = because the system is no more exhibiting 
a clustered phase. 

However, within a canonical ensemble it is more appro- 
priate to analyze M and P as a function of the temper- 
ature T — P^^. Both these quantities decrease with T 
and vanish when T > Tc ^ 0.55, their dependence from 
(Tc - r)-i is shown in Fig. |(b). Close to the critical 
temperature both M and P exhibit a jump indicating a 
first order transition between CP and HP. 

Thermodynamic potentials can be straighforwardly 
obtained from the partition function : 

hm f^logz)=^iM^ilM, 



U-- 



d{pF) dG{f3,UM,up) 



d(3 



d(3 



(18) 



where F — F{f3) is the Helmholtz free energy and U the 
specific (or internal) energy. In order to find the sta- 
ble configuration of the system we should look for the 
absolute minimum of F. The inverse temperature de- 
pendence of the minimum of F is reported in Fig. ^(a) . 
In the limit — > the free energy will have a minimum 
at the origin {uM,up) = (0,0) and the HP will be ob- 
served, while in the limit (3 —f oo the minimum will be 
located at < wp < um and the system will be in the 
CP. In the range f3 £ [1.81, 2.0], F displays two coexisting 
minima, that we denote, with straightforward notation, 
as Fhp and Fcp- However for /3 > /3c = 1.84 the CP is 
observed, because Fcp < Fhp, while for f3 < (3c the HP 
prevails since Fqp > Fhp- At (3 = (3c the two minima 
are equivalent {Fcp = Fhp), this is a further indication 
that the transition is of first order. 

Fig. ^(b) represents the temperature T as a function 
of the internal energy U. For small energies, U ^ T/2 



(i.e. there is a virial) and the particles are all trapped in 
a single cluster. For U > Uc ^ 2, the system is in HP 
and T is also increasing linearly with U. This indicates 
that the system behaves like a free particle gas. These 
two regimes correspond to the two integrable limits of 
model (0). 

In the intermediate energy range, the tendency of the 
system to collapse is balanced by the kinetic energy and 
T is no more proportional to U. For 1.6 < U < 2, the 
amount of kinetic energy is such that a significant frac- 
tion of particles escape from the cluster. The system 
is then characterized by the presence of two coexisting 
phases: the first is supported by the fraction of LEP 
that are trapped in the well of the self-consistent poten- 
tial and the second is due to the HEP that have large 
kinetic energy and thus behave almost freely. 



B. Microcanonical Results : Numerical Findings 

In order to obtain microcanonical results, standard 
molecular dynamics (MD) simulations have been per- 
formed for the present model (0) within a NVE ensem- 
ble [Q. The equations of motion are integrated start- 
ing from their formulation in terms of the single particle 
Hamiltonian hi (|^), where the dependence of the dynam- 
ics of each particle from mean-field variables is made ex- 
plicit. This approach turns out to be quite efficient from 
a computational point of view, because the CPU time 
increases only linearly with N (and not proportionally 
to N'^ as usual). This allows us to perform simulations 
of the model (|l|) for high numbers of particles (up to 
N ~ 10^) and for quite long integration time (up to 
t — 5 X 10^). As integration scheme we have adopted 
a recently developed 4*'* order simplectic algorithm [ p5[ : 
with an integration time step dt — 0.3 relative energy 
fluctuations AU/U remains smaller than 10^^ for any 
considered N and U. 

In the simulations, the particles are initially clustered 
in one single point and with a "water-bag" distribution 
for the velocities. Other initial conditions, for example, 
with particles uniformly distributed over the box and 
with Maxwellian velocity distribution have been tested. 
But as far as thermodynamical properties were concerned 
no significative differences have been observed. In or- 
der to avoid transient effects, the reported thermody- 
namical and dynamical quantities have been measured 
starting from initial " relaxed" states resulting from suffi- 
ciently long preliminary simulations. Moreover, the aver- 
ages have been typically performed over total integration 
times t ~ 10^. In the present section, the thermodynam- 
ical properties will be considered, while the dynamical 
ones will be the subject of the two following sections. 

We define the measured temperature T, following the 
equipartition theorem, through the time averaged spe- 
cific kinetic energy T =< K{t)/N >t. In Fig. |(b), 
the measured ?7-dependence of T is reported. The MD- 
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estimation of < M >t and < P >t are shown in Fig. ||(a). 
Both figures indicate that the numerical findings are in 
good agreement with the canonical predictions, apart in 
the interval 1.6 < J7 < Uc, where canonical estimate are 
not available. In this energy range (enlarged in the inset 
of Fig. |^(a)) we observe a decrease of T with U. This in- 
dicates that the specific heat will be negative for energies 
slightly below Uc- The existence of such a negative spe- 
cific heat regime in proximity of the " declustering" tran- 
sition for gravitational potential was predicted in 1971 by 
Hertel and Thirring These authors studied a simple 
classical cell model and noticed the non-equivalence of 
canonical and microcanonical ensembles in this region. 
These predictions have been successively confirmed by 
numerical findings for short ranged non-singular attrac- 
tive potential Pjq]. Negative specific heat regimes are 
forbidden within the Gibb's canonical ensemble, because 
physically unstable, and thus they are bridged by a con- 
stant temperature line (that in our case is T^. « 0.54). 
This picture is similar to the Maxwell construction for 
the van der Waals isotherms in the liquid-vapor coexist- 
ing region. This prohibition does not hold in the micro- 
canonical ensemble and the physical implication of this is 
fundamental in astrophysics where negative specific heat 
regimes have been studied for several decades [Q. 

Within the microcanonical context, the phenomenon 
of negative specific heat can be understood following a 
heuristic argument [|l[. As shown in Fig. ^(a) approach- 
ing the transition the values of M and P decay very fast. 
This suggests that a very limited increase of U yields a 
significant reduction of the number of LEP and thus a 
strong increase of the potential energy. Total energy be- 
ing conserved, this excess of potential energy has to be 
compensated by a loss of kinetic energy, as a result the 
system will become cooler Q]. 

From our numerical investigations it turns out that for 
any U the average thermodynamical quantities (namely, 
T,M and P) are independent from N (for 100 < TV < 
10, 000), apart at the critical energy. At U = Uc, a clear 
decrease of M is observed from a value 0.198 at iV = 200 
to 0.063 at TV = 10, 000 indicating that the system is ap- 
proaching the HP (a similar behaviour is observed for P). 
This is confirmed from the fact that also the temperature 
tends to the corresponding HP value {T — U — 3/2) for 
increasing TV: T varies from 0.55 at TV = 200 to 0.51 at 
TV = 10,000. This effect is probably due to the strong 
fluctuations in the single particle potential Vi due to finite 
TV-effects, that become dramatic at the transition. For 
small TV- values a clustered situation seems to be favorite, 
while in the mean-field limit the HP is finally recovered. 



IV. DYNAMICAL PROPERTIES 



A. Dynamical transition and transport properties 

Once examined the microcanonical thermodynamical 
results, we now concentrate on the description of the 
transport properties of (|l|). In the CP each particle of 
our system moves in a single particle potential Vi = Vi {t) 
(^) that instantaneously is similar to the so-called egg- 
crate potential (see Fig. 0) ||2|,|l^,||-f| . The study 
of single particle motion in egg-crate fixed potential was 
originally motivated in the context of adatoms diffusion 
over a rigid surface |36| . For a single particle moving 
in such fixed potential the self dynamics is known to be 
anomalously diffusive when the particle is channeling, i.e. 
when its energy lies between Vm and Vs | |2^ , p^ , p3[ -p5[ . 
This anomalous behaviour is due to the competition of 
laminar and localized phases. Indeed, channeling parti- 
cles intermittently show an almost ballistic motion along 
the channels of the potential interrupted by localized se- 
quences, where the particles bounces back and forth on 
the maxima of the potential. 

Usually anomalous diffusion has been studied in sys- 
tems with very few degrees of freedom |l^, 0. Only few 
att emps h ave been made to consider TV-body dynamics 
||l3| , p7 28 1 . One of the reasons of this is that the theoret- 
ical analysis of dynamical properties of high dimensional 
Hamiltonian systems is particularly complicated and not 
fully understood. Moreover, accurate simulations of such 
systems for TV >> 1 and over long time intervals, neces- 
sary to understand the asymptotic diffusive regime, are 
quite difficult to perform. 

In this section we first focus on the numerical descrip- 
tion of the transport properties of (Q) and we will try 
to give some indications relative to the basic dynamical 
mechanisms governing single particle transport. To this 
aim, we consider the time dependence of the mean square 
displacement (MSQD), that usually reads as 



(19) 



where the average < . > is performed over different time 
origins and over all the particles of the system. The trans- 
port is said to be anomalous when a ^ I: namely, it is 
subdiffuse if < a < 1, superdiffusive if 1 < a < 2 and 
ballistic for a = 2 1 17 2^Jl^ . The usual Einstein diffusion 
law corresponds to a = 1 and in 2D can be written as 
< r^(t) >= ADt, where D is the self-diffusion coefficient. 
We consider in this paper, D and a as the basic relevant 
observables for the description of transport. 

More refined diagnostics might be considered, for ex- 
ample the probability distribution of the time inter- 
vals whithin which the trajectory of the particle re- 
mains trapped. The latter was explicitely computed 
both for simple maps and for a single particle moving 
in a fixed egg-crate potential landscape and shown to 
exhibit a power-law deca y |]l7|,|l8| , p2 37 3^ responsible of 
the anomalous diffusion p2| , |23{ . In the limit TV ^ oo 
we expect to find similar indications also for model (|^) 
once the lEP's dynamics is considered. But for finite TV, 
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the energy of each particle will not be constant in time. 
Therefore a particle initially of type lEP can become a 
LEP or a HEP, giving rise to a much more complicated 
dynamical behaviour, that cannot be simply described 
in terms of the localization time distribution. However, 
if one considers an initially lEP and follows its trajec- 
tory for some time it displays features quite similar to 
the so-called Levy walks 18 1 (see Fig. g). This kind 
of trajectories has been usually identified when anoma- 
lous diffusion occurs. Therefore, we believe that also for 
model (|I]) anomalous diffusion will be observable. 

As we already mentioned, finite iV-effects play a de- 
terminant role in the dynamics. Due to self-consistency, 
they are responsible for the fluctuations in time of the 
mean-field quantities Mi^2 and Pi, 2- The potential ex- 
perienced by each particle thus fluctuates in time and 
particles having an energy close to Vg have the possi- 
bility to be trapped in the potential well as well as to 
escape from it. This implies that the localization phe- 
nomena illustrated in Fig. ^ are not only due to bounces 
of the particle on the maxima of the potential, but also to 
trapping in the potential well due to separatrix crossing. 

We have argued from simple considerations that dif- 
fusion should be anomalous in CP. This is indeed the 
case, as confirmed from direct evaluation of the MSQD 
in a quite extended energy interval. An example of this 
is reported in Fig. | (a) for J7 = 1.1 and TV = 4, 000. The 
diffusion is anomalous for times smaller than a crossover 
time T beyond which the Einstein's diffusion law is re- 
covered < r^{t) >= 4:Dt. A similar behaviour for the 
MSQD was previously observed for a system of N sim- 
plectic (globally and locally) coupled maps 28|, but 
with a subdiffusive (i.e. with a < 1) short time dynam- 
ics. 

The direct study of the velocity autocorrelation func- 
tion (VACF) Z{t) confirms the general features seen for 
the MSQD in the CP: namely, on times t < ty the VACF 
is characterized by a long-time tail that decays as t°'~^. 
This power-law decay is fully consistent with the cor- 
responding one observed for the MSQD For times 
longer than , the VACF decreases exponentially as usu- 
ally expected for brownian motion (see Fig. ^ (b)). It is 
therefore reasonable to expect that t (x t^. 

The energy dependence of the a-values is illustrated 
in Fig. 0. It shows up clearly that the thermodynamical 
phase transition from CP to HP is associated to a dy- 
namical transition from superdiffusion (with 1.3 < a < 
1.9 for OA<U < 2.0) to ballistic motion (with a ~ 2 for 
U > Uc)- In the CP regime, we observe an increase of a 
from 1.3 ±0.1 to 1.9 ±0.1, that is due to the modification 
of the shape of the single particle potential. 

It is reasonable to expect that this enhanced diffusive 
behaviour with an energy dependent exponent is linked 
to the fraction of channeling particles. The number of 
such particles is related to the depth of the potential 
well (namely, to Vm — Vm) and to the energy width of 
the channels Vm — Vs. At very low energy U < 0.3, all the 
particles are essentially trapped in the potential well and 



none is channeling: no diffusion is observed in this case. 
For increasing energy a fraction of particles (due to the 
decrease of Vm ~ Vm) will escape from the well and some 
of them get enough energy to move along the channels: 
anomalous diffusion is then evidenced. The increase in 
the value of the exponent a is due to the fact that also 
the channel width Vm — Vs grows with U (see Fig. |^). 
However, for U approaching Uc the number of untrapped 
particles increases noticeably, but now the channel width 
vanishes abrubtly this implies that a significant fraction 
of particles will move freely (with energy > Vm)- These 
mechanisms lead naturally to ballistic motion for U > 
Uc, where the potential Vi is now almost constant apart 
fluctuations of order 0{1/\/N)- 

The fact that in the asymptotic limit [t — > oo) normal 
diffusion is recovered constitutes a typical signature of a 
noisy dynamics 1 29 3^,^ . In order to avoid artifacts due 
to numerical noise, in the iroplementation of the integra- 
tion scheme for the model (|l) we took care to maximixe 
the numerical precision. Therefore, the transition from 
anomalous to asymptotic ordinary diffusion is attributed 
to a " deterministic source of noise" , that is intrinsic of 
our system and due to finite size effects, as we'll show in 
the next section. 



B. Finite effects 

In this section, we consider finite N effects in order 
to understand the asymptotic time dependence of the 
MSQD outlined in the previous section. As a matter of 
fact, for finite A^, due to the crossing from anomalous 
to normal diffusion on long times, a standard diffusion 
coefficient D can be always defined in the limit i — > oo. 
Therefore we will focus our analysis on the A^-dependence 
of D and of the crossover time r. 

A deterministic anomalously diffusing dynamical sys- 
tems submitted to weak environmental white noise shows 
a transition from anomalous diffusion to standard diffu- 
sion on sufficiently long time scales The crossover 
time T for this transition is explicitely computed and 
turns out to increase as an inverse power of the noise am- 
plitude, when short time behaviour is superdiffusive. In 
this context, D is shown to be directly proportional to a 
certain power of r | |30[ |. 

Starting from this knowledge, we can understand the 
nature of the mechanism that generates asymptotic stan- 
dard diffusion in model (^. In order to give an un- 
ambigous definition of t, we consider the local slope of 
ln(< r'^{t) >) as a function of ln{t)- The crossover time 
T is determined when this slope becomes smaller then a 
threshold value /i. The diffusion coefficient D can be 
related to r, assuming that t oc ty (as already men- 
tioned). By the definition of the diffusion coefficient we 
have Z? cx Z{t)dt. Assuming that ty is sufficiently 
long, the following relationship is then straightforwardly 
found 
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D cx r"-\ (20) 

This result is in agreement with theoretical predictions 
and is sucessfuUy confirmed by numerical simulations of 
noisy maps pO[ . 

In Fig. ||we report the N dependence of r for /i = 1.1 
and for two energy values, namely U = 1.48 and U — 2.0. 
In both cases, we find that t oc N. Moreover, this depen- 
dence is not related to the chosen value for the threshold 
fi. We have indeed verified that for fi = 1.2 no quali- 
tative difference could be detected. The interpretation 
of the N dependence of the crossover time is straighfor- 
ward if we consider finite N effects as a source of noise 
in our model and we can directly deduce that the ampli- 
tude of the noise is typically of order 1 / y/N. This last 
assumption is justified by the fact that the microscopic 
dynamics of the particles generate stochastic fluctuations 
0{1/^/N) in the values of M and P. These fluctuations 
become weaker for increasing N and thus naturally yield 
an increasing value of r. 

The idea of a dynamical noise source due to finite TV is 
also confirmed by the analysis of D, that turns out also to 
be an increasing function of N. Some data are reported 
in Fig. I for J7 = 1.48 and U = 2.00. In particular, 
considering systems with 100 < < 10, 000 we find D cx 
iVT with 7 = 0.7±0.1 and l.OiO.l for U = 1.48 and 2.00, 
respectively. These numerical evidences suggest that 

D cx V-'< (21) 

where V is the variance of the white noise applied on 
the system. For model (|l|), we have V cx 1/7V. Being r 
proportional to N, from Eq. ( |2C| ) we expect that 

7 = a-l . (22) 

However, we find also a dependence of the a- values on 
N, not neghgible at least for iV < 3, 000. But for N > 
3, 000 a saturation to asymptotic values is finally achieved 
(as clearly shown in Fig. |ll| ). Assuming for a their 
asymptotic values, we obtain from the relation (^2|) the 
following 7- values: ~ 0.64 for U = 1.48 and ~ 0.9 for 
U = 2.00. In view of the finite N limitations, we consider 
these results as consistent with the numerical results. A 
relation analogous to (|2^) was previously found in | |4C| ] 
for the eddy diffusivity associated to a three dimensional 
noisy velocity field Moreover, the authors of Ref. 

pof have shown that Eq. (|2^) holds also for their model. 
Due to the complete different nature of the two system, 
we expect that Eqs. ( ^l| ) and (|2^) should have some more 
general field of applicability. 

In the context of high dimensional Hamiltonian sys- 
tems, previous results showed that the crossover time 
to normal diffusion is inversely proportional to the dif- 
fusion coefficient when short time behaviour is subdif- 
fusivc 1^,^. The crossover to standard diffusion is 
then interpreted as a consequence of the destruction of 
the self-similar structure of the stability islands in phase 
space. We believe that the asymptotic normal diffusive 



behaviour that we observe has the same origin as in [£7| . 
But let us discuss more in detail the dynamical mecha- 
nisms that are present in our case. 

For finite N, due to the self-consistent character of 
model (|l|) , the single particle potential fluctuates in time 
with typical amplitude 0{1/^/N). The values of the sad- 
dle points thus also fluctuate in time and naturally gen- 
erate a time pulsating separatrix sweeping a phase space 
domain of width 0{1/Vn). Hence, a particle with en- 
ergy close to Vs can cross the separatrix and stochasti- 
cally experience trapped and channeling motions. We 
showed that transport is anomalous in CP and relies 
on the channeling particles that exhibit localized motion 
when they bounce back and forth on the maxima of the 
potential. But, for finite N, the pulsations of the separa- 
trix induce a second localization mechanism, that tends 
to trap the trajectories of lEP into the potential well. 

To better illustrate the different dynamical behaviours 
occuring in the present system, let us focus on the dy- 
namics of three typical particles, namely an initially HEP, 
lEP and LEP. We register for each of them the time 
evolution of (a) their orbit in the (x, y) plane, (b) their 
coordinate x{t) and (c) the single particle energy hi as 
defined in (||). The results are shown in Fig. |l^ where the 
first column corresponds to HEP, the second to lEP and 
the third to LEP. For each particle (a), (b) and (c) are 
plotted in the first, second and third row, respectively. 
Visual inspection of HEP (a) and lEP(a) shows a quiet 
similar behaviour and an enlargement of these trajecto- 
ries would yield a picture similar to the one already re- 
ported in the inset of Fig. ||(a). However, Fig. [Tc|-IEP(b) 
indicates that the initially lEP experiences long time lo- 
calized sequences that are not present for the HEP. Fig. 
|l^IEP(c), showing the time dependence of the energy 
of lEP, clearly indicates that in the long lived localized 
regime corresponding to t € [10000,17000], the particle 
is trapped since its energy hi goes below the average en- 
ergy of the separatrix < Vs > (indicated in this figure by 
a full line). 

The localization within the potential well, due to sep- 
aratrix crossing, interrupt the sequence of correlated 
flights and localizations along the channels, and naturally 
inhibits superdiffusion. This inhibiting effect of separa- 
trix crossing is due to finite N and its consequences on 
the numerical value of the exponent a is illustrated in 
Fig. |ll| for two values of U. For given U, this figure 
shows an increase of the order of 10% in a when rais- 
ing N from 10'^ to 10'*. The exponent a then saturates 
to an almost constant value indicating that the effect of 
trapping and untrapping by the fiuctuating separatrix de- 
creases and finally becomes negligible for growing N . In 
other words, as N grows, the phase space volume swept 
by the pulsating separatrix shrinks and finally vanishes 
in the mean-field limit N ~^ 00. As a final remark, we 
should also notice that the channeling particles can be 
" decorrelated" not only by the trapping mechanism but 
also by an escaping mechanism to higher energies (for 
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hi > Vm ) ■ But from numerical evidence we can conclude 
that this other mechanism is less relevant, at least on the 
time scale of our simulations. 

Finite TV-effects are responsible, besides of the fluctu- 
ations of the single particle potential, also of a second 
mechanism generating "dynamical noise". In order to 
explain it more in detail let us consider the LEP particle 
dynamics: a typical LEP-orbit is reported in Fig. 
One observes essentially two main aspects of its dynam- 
ics: (i) the extremely slow motion of LEP with respect 
to lEP and HEP (cf . the scales) ; (ii) a drifting motion of 
LEP that resembles a brownian motion. During all the 
simulation time (t = 30, 000) the corresponding energy, 
plotted in Fig. ra-LEP(c), remains close to the average 
value of the potential minimum < Vm >~ 0.654. This 
indicates that the particle remains trapped during the 
entire simulation. Being for U — 1.0 the crossover time 
10'' < r < 10^, the previous result suggest that such a 
trapped particle will eventually escape from the well on 
time scales much longer than r. We moreover verified 
that this observation remains true for almost any parti- 
cle with initial energy smaller that O.bVg « 1.85. Which 
for U — 1.00 represent approximatively 75% of all the 
particles. 

The above observations confirm the fact that anoma- 
lous diffusion is only due to the small fraction of parti- 
cles that evolve inside the channels of the potential. The 
trapped particles, indeed, due to their extremely slow 
motion, contribute negligibly to the average MSQD. But 
the slow drift of the clustered particles, corresponding to 
a drift of the potential seen by a single particle, inhibits 
fiights on very long times. This because a channeling 
particle will not have a free horizon in front of itself at 
any time, as it happens for a particle moving in a fixed 
potential frame. 

The origin of the collective drifting motion of the clus- 
tered particles can be understood in the following way: 
let us assume that at time t = all the particles are 
trapped in the potential well and that the total momen- 
tum is zero and that at a later time t > one particle 
escapes from the cluster. This particle will carry out of 
the cluster a non-vanishing average momentum p, but 
as the total momentum should be conserved, the clus- 
ter will then start to move with an average momentum 
—p/{N— 1). Therefore, in the mean-field limit this effect 
will disappear and consistently we will have anomalous 
diffusion at any time . 

The two mechanism outlined above, will have a cumu- 
lative effect on the phase space topology of our system. 
Due to the analogy with the egg-crate potential, in the 
mean-field limit the phase space of our model will ex- 
hibit a hierarchy of nested self-similar stability islands 



lands disappear up to a typical size that decreases for 
growing N. Being the self-similarity no more complete 
in phase space, normal diffusion should be recovered be- 
yond a crossover time that grows for increasing N. 



V. LYAPUNOV ANALYSIS 

In order to complete the description of model (j^), we 
investigate in this section a fundamental indicator to 
characterizes the dynamics of Hamiltonian models: the 
maximal Lyapunov exponent A. Our analsys rely on nu- 
merical estimation of A, performed considering the evo- 
lution in the tangent space of the model and applying a 
standard technique introduced in Ref. . 

Our model is integrable in the limit of low and high en- 
ergy, therefore A — > for [/ — > and J7 — > oo. In between 
these two extrema we expect that a finite Lyapunov ex- 
ponent will be observed similarly to what was recently 
found for ID self-gravitating models |43-45|. Our data 
are reported in Fig. 12 for three different types of initial 
conditions : 
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^ ^ ^ , 25[1- The monotonous increase of the value of 

a with N reported in Fig. |ll| intuitively suggests a simi- 
lar continous picture when considering the modifications 
occuring in the phase space topology under continous 
variation of N. In particular, for finite N we expect that 
the phase space structures associated to the smallest is- 



• (A) the particles are initially clustered and the ve- 
locity distribution is Maxwellian; 

• (B) the velocity distribution is again Maxwellian 
but with a thermal velocity conciding with its 
canonical prediction and the particles are organized 
in a single cluste in such a way that also the M- and 
P- values coincide with their canonical prediction; 

• (C) the particles are initially clustered with a 
water-bag velocity distribution. 

The initial condition (C) is the one commonly used 
through the present paper, in particular for the study of 
transport in the system. As can be seen fom Fig. ^ A 
grows for increasing [/ up to a maximum value and then 
decreases. Such maximum (at least for N — 200) is lo- 
cated at an energy U ~ 1.3 — 1.4 < Uc- For U > Uc, we 
observe a power law decrease of A with N. In particu- 
lar, for U = 3.0 we found a power law exponent ~ 0.31 
(see Fig. |l^, in good agreement with recent theoretical 
results obtained via a random matrices approach in Ref. 
iQ and with a Riemanian geometrical technique p5[ . 
In both these studies an exponent 1/3 has been found 
considering a model similar to (|I]) in ID 0. 

In the low energy limit (for U < 0.01) a power law 
increase of the type U^^'^ is clearly observable for all the 
three types of initial conditions. The data are reported 
in the inset of Fig. A similar behaviour was found for 
a ID mean-field model [Q. This indicates that for fully 
coupled Hamiltonian systems this property holds in gen- 
eral, independently of the space dimensionality. In par- 
ticular, for [/ ^0 the particles are all clustered, therefore 
we expect that the scaling A oc C/^/^ should be related to 
a "collective" chaotic mechanism. Work is in progress 
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in order to derive a theoretical explanation of such be- 
haviour p6| . 

Let us now try to understand which are the mecha- 
nisms underlying the observed behaviours of A. Once 
the energy U is fixed for all the three considered ini- 
tial conditions, after a reasonable transient, we obtain 
exactly the same value for the average temperature T 
and magnetizations M and P. However, at low ener- 
gies {U < 0.8) the measured average A depend heavily 
on the initial conditions. This clearly indicates the co- 
existence of several equivalent state, that can be consid- 
ered as equilibrated within the examined time interval. 
Usually we have averaged the maximal Lyapunov over 
a time 1,000,000 < t < 9,000,000 after an equilibra- 
tion time ranging from t = 500,000 to t = 20,000,000 
(this last value has been used in particular for extremely 
low energies). Obviously, also if the considered time scale 
are considerably long we cannot exclude that these states 
are metastable. It should be noticed that this kind of be- 
haviour is unexpected in TV-body Hamiltonian systems, 
because it is commonly believed that for sufficiently high 
N Arnold diffusion takes place and each orbit is allowed 
to visit the complete phase space. But our data instead 
indicate that some "barrier" in the phase space still sur- 
vive even for N = 200. The origin of this lack of ergodic- 
ity is related to the long range nature of the forces that 
induces a persistent memory of the initial conditions, as 
previously noticed by Prigogine and Severne for grav- 
itational plasmas. Recently, some numerical evidence of 
non ergodicity has been reported also for a ID mass sheet 
model ]4^ ]. 

It is clear from Fig. [l^that, in the interval U G [0, 0.8], 
for initial conditions of type (B) A remains always smaller 
than the corresponding exponents obtained with initial 
conditions (A) and (C). The maximal differences are ob- 
served in the energy range 0.2 < U < 0.8, where particles 
begin to escape from the cluster (this for initial conditions 
(A) and (C)). Above U ~ 0.9 the same Lyapunov is ob- 
tained for all type of initial conditions. A typical feature 
of the initial conditions (B), for U < 0.8 is that all the 
particles are trapped in the potential well. Instead when 
one or more particle escape from the cluster {U > 0.9) 
also with this initial condition the usual A is obtained. 
We believe that two chaotic mechanisms are present in 
the system: one felt from the particles moving in the 
minimum of the potential and one from particles visiting 
a region near to the separatrix. This second mechanism 
is known and is related to a chaotic belt present around 
the separatrix pgf . 

The first one should be instead due to the temporal 
erratic oscillations of the minimum of the potential well. 
In order to clearly identify such mechanisms we have fol- 
lowed the trajectory of a system initially prepared with 
condition (B) for an energy U = 0.87. On short time 
all the particles are trapped and we measure an average 
value A ~ 0.13. At a later time one particle escapes from 
the cluster and A shows a jump to a value that is al- 
most double (see Fig. |lj). In Fig. |lj the magnetization 



M and the kinetic energy K are also reported. When 
the particle escapes from the cluster M shows a clear 
decrease as well as K. This last effect is due to the fact 
that the potential energy V is minimal when all the parti- 
cles are trapped, therefore if one escapes V increases and 
due to energy conservation K decreases. This is the phe- 
nomenon at the basis of the negative specific heat effect. 
From the above arguments we can identify a strong chaos 
felt from the particles approaching the separatrix and a 
minimal chaos associated to orbits trapped in the min- 
imum of the potential. The presence of these 2 chaotic 
mechanisms together with the non-ergodicity of the sys- 
tem explains the strong dependence of the A- values from 
the starting conditions. 

As a final point, we would like to notice that in the 
CP for low energy density the Lyapunov exponent aver- 
aged over short times exhibits an intermittent behaviour, 
when starting conditions (A) or (C) are considered. In 
particular in Fig. |l^ the instantaneous A is shown for 
condition of type (C) together with its running averaged 
value and the corresponding averaged value for condition 

(B) . The intermittent behaviour can be explained as due 
to trapping and untrapping of the lEP's. As a matter 
of fact, for initial conditions (B) and energies U < 0.8 
no diffusion at all is observed (the MSQD saturates to a 
constant value for long times), while for initial condition 

(C) (or (A)) a super-diffusive motion is observed for ener- 
gies higher than U = 0.3. This confirms that channeling 
particles (i.e. lEPs) are affected by a stronger chaotic 
mechanism than the trapped ones (namely, the LEPs). 

VI. CONCLUSIONS 

In the present article we have analyzed the equilibrium 
and dynamical properties of a 2D A^-body self-gravitating 
Hamiltonian system. Our main result is the occurence of 
a thermodynamical phase transition associated to a dy- 
namical transition from anomalous to ballistic transport. 

Firstly, the statistical equilibrium description in the 
canonical ensemble has been reported. Revealing a first 
order phase transition between a collapsed phase, char- 
acterized by the presence of a single cluster of particles, 
and an homogeneous phase where particles are uniformly 
distributed. Close to the transition, the presence of two 
equivalent minima in the Helmholtz free energy clearly 
indicate a regime where the two phases (CP and HP) 
coexist. Within the microcanonical ensemble, for ener- 
gies close to Uc MD results reveal a negative specific heat 
regime. Referring to previous work ||l|, we interpret this 
phenomenon as a signature of non-equivalence between 
canonical and microcanonical statistical ensembles. 

Secondly, we have shown that in a broad interval of 
energies 0.3 < U < Uc the transport in this system is 
anomalous. In the sense that the particles show a super- 
diffusive motion below a cross-over time and normal dif- 
fusion for longer times. Usually, normal diffusion is ex- 
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pected to occur in iV-body systems due to the absence 
of long time correlations. Instead in the present case the 
non markovian nature of the process is clearly evidenced 
by the occurence of Levy walks and of power law decay 
for the VACF. In the CP, from simple mean- field consid- 
erations one can conclude that each particle evolves in a 
2— D time dependent egg-crate like potential. Channeling 
particles are at the origin of anomalous diffusion on short 
time scales. Since their motion is made up of flights along 
the channels interrupted by localizations due to bouncing 
of the particles on the maxima of the potential. We have 
also shown that finite N effects generate a pulsating sep- 
aratrix and give rise to a second localization mechanism 
for the particle motion due to the trapping of the orbits 
in the potential well. This phenomenon together with 
the slow drift of the potential well is responsible for the 
suppression of anomalous diffusion on long time scales. 

These two finite A^-effects can be interpreted as a white 
noise source affecting the single particle dynamics. As a 
matter of fact these effects become weaker for increas- 
ing N, as confirmed by the linear dependendence of the 
crossover time on A^. The asymptotic dynamics of the 
model then depends on the order the two limits A — s- oo 
and t — > oo are taken. If we first take the limit N oo 
and then t oo diffusion remains anomalous for all time. 
Ohterwise standard diffusion is recovered for sufRciently 
long time. 

The Lyapunov analysis shows clearly that even for 
N = 200 the phase space is not fully accessible and the 
visited portion of the phase space depends on the initial 
conditions. This is quite unexpected in high dimensional 
Hamiltonian systems. Moreover, two different chaotic 
mechanisms have been identified: a minimal chaos as- 
sociated to the clustered particles and a strong chaos af- 
fecting the lEP particles. Anomalously diffusing particles 
have therefore an intermittent chaotic behaviour. 

We believe that the properties displayed by the model 
(|l|) are quite general for gravitational systems. It is also 
our claim that the inclusion of higher order Fourier har- 
monics in the potential will not qualitatively affect the 
results here presented. This is indeed confirmed by a re- 
cent study by Abdalla and Reza Rahimi Tabar for 
the 2D full logarithmic Newtonian potential Voo- Also 
in that case a transition at finite T from a clustered to 
a homogeneous phase is still present. This constitutes 
in our opinion an interesting basis for the generalization 
of the present results to models that mimic more realis- 
tically the gravitational interactions (e.g. systems with 
true Newtonian potentials in 3D). 

Due to the strict analogy of the model here studied 
with point vortices model for 2D turbulence |^l| we ex- 
pect that anomalous diffusion should also be observable 
for such models. Preliminary indications | pT[ | seem to 
confirm our claim and ask for more accurate investiga- 
tions in such direction. 
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Fig. 1 - Torcini et al. 



FIG. 1. Instantaneous single particle potential V{x,y) 
together with its contour-plot for U = 1.00 and 
{x,y) e [-57r/2,57r/2]^ Mi = Ma = 0.8249, 
Pi=P2 = 0.6929, (t>z=G = il}^ 
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FIG. 2. P dependence of um (full line) and up (dashed 
line). For /3 € [1.81,2] there are two coexisting minima: one 
at the origin {um,up) = (0,0) and a second one for which 

< Mp < UM. 
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FIG. 3. In (a), the canonical estimation of M and P (solid 
curve) together with their time averaged values obtained from 
microcanonical simulations (symbols) are shown. The mea- 
surements have been obtained with N = A, 000 (a part a few 
point with A'^ — 10*) and averaged over a total integration 
time ranging t = 1.2 x 10*' to t — 2.4 x 10^ with a time step 
dt = 0.3. In (b) arc reported the canonical estimation for M 
and P as a function of {Tc — T)~^ in a log- linear plot. Neither 
M nor P have usual critical exponent. 



FIG. 4. In (a) we plot the minimum value of the free 
energy function F as a function of /3. In the coexisting range 
F has two minima which take the same value at critical inverse 
temperature /3c « 1.84 (see inset), (b) Temperature as a 
function of the specific energy U . The theoretical predictions 
obtained within the canonical ensemble are indicated by lines. 
The circles correspond to microcanonical numerical results. 
The inset is an enlargement of the transition region: the full 
curves (resp. dashed) refer to the absolute (resp. relative) 
minimum of F{T). 
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Oe+00 le+04 2e+04 3e+04 log (t) 
t 

FIG. 5. Typical orbit of an initially lEP in the {x,y) plane FIG. 6. Time dependence of the mean square displace- 

(a) for U = 1.00 and N = 4,000. In this representation, the ment < r^(t) > in a log- log representation for U — 1.1 and 

2 — D torus on which the dynamics takes place is unfolded A'^ = 4, 000. The numerical results correspond to the sohd 

and represented as an infinite plane constituted of an infinite line. The segments are the estimated slopes of < T'^(t) > 

number of juxtaposed elementary periodic cell of size 1x1. In when t < t and t > t. In (b) the logarithm of the VACF is 

the inset an enlargement of the trajectory in the indicated box displayed as a function of log(t) for U = 1.74, in such a case 

is reported. In (b) we plot the corresponding time evolution the estimated a-value is 1.59. The reported slope is 2 — a = 

of coordinate x{t) and the enlargement corresponding to the 0.41 
sequence in the inset of (a). 
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FIG. 7. Coefficient a (triangles), time averaged depth of 
tfie potential well Vm — Vm (diamonds) and of the channels 
Vai — Va (circles) of the single particle potential as a function 
of U. The full line corresponds to the canonical prediction. 
The measurements of a have been obtained with A'^ = 4, 000 
(a part a few point with A'^ = 10^) and averaged over a total 
integration time ranging t — 1.2 x 10^ to t — 2.4 x 10® with a 
time step dt = 0.3. 
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FIG. 9. Logarithm of the diffusion coefficient D as a func- 
tion of log(iV) for U = 1.48 (diamonds) and U = 2.00 (circles). 
Straight lines correspond to linear fitting to the numerical 
data, the slope are 0.7 ± 0.1 and 1.0 ± 0.1 for U = 1.48 and 
U = 2.00, respectively. 
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FIG. 8. Crossover time r as a function of TV for U — 1.48 
(circles) and [7 = 2 (diamonds). The values of r have been 
estimated when the local logarithmic slope of the MSQD de- 
creases below a threshold value /i — 1.1. The solid line is 
the best linear fit to the data. Its slope is 0.95 ± 0.08 (resp. 
0.96 ± 0.06) for U = 1.48 (resp. U = 2.00). 
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FIG. 10. Typical orbits of HEP, lEP and LEP (first row) 
together with the time dependence of coordinate x{t) (second 
row) and the time dependent energy hi (third row) of each 
orbit. The first column refers to HEP, the second to lEP 
and the third to LEP. Note the difference of scale in graphs 
of the first and second line. N = 4000 and U = 1.00 and 
t < 30000. In figure lEP(c), the straight line is the average 
energy < Vs >~ 3.694 of the scparatrix and the dashed one 
the average energy < Vm >~ 3.958 of the maximal value of 
the potential. 



FIG. 12. Maximal Lyapunov exponent as a function of U 
in a log- linear plot, for three different types of initial condi- 
tions (for details see the text). The dashed line indicates the 
critical energy Uc- In the inset, log(A) versus log(f/) is re- 
ported at low energy. A scaling A cx U^^^ is evident in the 
low energy limit for all the three initial conditions. All the re- 
ported data corresponds to iV = 200 and to integration times 
1 X 10® < t < 9 X 10®. 
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FIG. 11. Exponent a as a function of N for U — 1.48 
(diamonds) and U — 2.00 (circles). Dashed lines correspond 
to exponential fitting. The data have been obtained for a 
total integration time t = 1.2 x 10®. 
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FIG. 13. log(A) versus log(Af) for = 3.0 > Uc- 
straight line is a linear fit to the data with slope 0.31. 
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FIG. 14. Time evolution of the maximal Lyapunov expo- 
nent, of the magnetization M and of the kinetic energy K are 
reported for an initial condition of type (B) for \J = 0.87 and 
N = 200. 
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FIG. 15. The maximal Lyapunov exponent averaged over 
short times {t = 300) is reported for an initial condition of 
type (C) (solid line) together with the corresponding running 
average (long-dashed line). The lower dot-dashed line corre- 
sponds to the running average of A for an initial condition 
(C). The energy for both cases is U = 0.5 and N = 200. 



17 



